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Introduction 

The  purpose  of  this  paper  is  to  develop  and  to  illustrate  by  a  number 
of  examples  some  interesting  and  useful  extensions  of  a  procedure,  gener- 
ally known  as  Dunkerley's  method,  for  estimating  the  gravest  frequency  of 
an  undamped  linear  mechanical  vibrating  system.  The  theory  is  developed 
in  the  next  section  hereof.  The  equivalent  of  equation  15  was  given  in 
the  year  189^  by  S.  Dunkerley  (3)  who  obtained  the  result  empirically 
based  on  calculations  performed  while  investigating  the  vibrations  of 

shaft  and  disk  systems;  he  certainly  did  nofe  recognize  the  generality  or 
the  theoretical  basis  of  his  formula.  Temple  and  Bickley  (12)  discussed 

the  procedure  in  1933  indicating  Its  applicability  to  both  lumped  and 
continuous  mass  distributions.  However  Bickley  and  Talbot  (1)  In  a 
later  (1961)  textbook  on  vibrations  do  not  mention  the  method.  South- 
well (11)  treats  the  method  in  his  well  known  treatise  (1936,  19^1)  without 
citing  Dunkerely.  This  may  be  the  reason  that  the  method  is  sometimes 
called  the  Dunkerley-Southwell  method.  However,  although  Southwell  may 
have  arrived  at  the  result  independently  (he  did  indeed  discover  another 
method  of  determining  lower  bounds;  cf.  Lamb  and  Southwell  (6)),  his  book 
does  cite  Temple  and  Bickley  (12)  among  the  general  references. 
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A  formulation  in  terms  of  integral  equations  and  their  eigenvalues, 
iterated  kernels,  etc.,  which  is  appropriate  for  continuous  mass  distrib- 
utions, is  a  logical  extrapolation,  by  analogy,  from  the  formulation  in 
terms  of  matrices,  their  eigenvalues,  powers  of  matrices,  etc.,  which  is 
appropriate  for  systems  having  a  finite  number  of  degrees  of  freedom,  so 
that,  in  a  sense,  all  the  work  reported  herein  stems  from  Dunkerley  (3) 
and  Temple  and  Bickley  (12).  Indeed  it  was  via  this  path  that  the  writer 
was  led  to  the  results  reported  here.  However,  the  relationship  with  the 
theory  of  integral  equations,  a  thoroughly  developed  and  explored  discipline, 
is  so  very  close  that  it  seemed  likely  that  the  developments  at  which  the 
writer  has  arrived  had  been  anticipated,  in  purely  mathematical  context,  by 
an  earlier  writer.  Indeed  this  is  the  case.  Mikhlin  (8)  definitely  attri- 
butes the  ideas  to  Mikhlin  (7).  In  neither  of  these  references,  however, 
does  Mikhlin  refer  to  Dunkerley  (3),  Temple  and  Bickley  (12),  or  South- 
well (11).  It  seems  clear  that  Mikhlin  arrived  at  his  results  without 
being  aware  of  Dunkerley' s  formula  or  of  its  relation  to  his  own  work. 

Accordingly,  the  writer  believes  that  it  is  appropriate  to  call  the  ex- 
tended procedure  by  the  name  Dunkerley-Mikhlin.    The  present  paper  may 
be  considered  as  a  brief  exposition  of  this  method,  and  some  extensions, 
and  of  its  application  to  a  variety  of  problems  of  engineering  vibration. 

Derivations 

We  consider  an  undamped  linear  vibrating  system  characterized  by 
N*N  symmetrical  matrices  K  (stiffness)  and  M  (mass).  Suppose  that  there 
are  p  rigid  body  modes  which  are  known  (by  inspection).  (Thus,  if  p  >  0, 
K  has  no  inverse.)  There  exists  a  modal  matrix  U,  not  necessarily  unique, 
and  a  (diagonal)  spectral  matrix  ft2  in  which  frequency-squares  are  arranged 
in  non-descending  order,  such  that 
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KU  =  MUft2;     UTMU  =  I  =  N*N  unit  matrix  (la,b) 

Now  consider  any  other  system  which  is  similar  to  the  original  system 
except  that  p  additional  constraints  have  been  incorporated  so  as  to  eliminate 
the  rigid  body  modes.     Let  the  flexibility  of  this  system  be  C. 

Consider  the  kth  mode  (lop)  of  the  original  system.     There  are  constants 
a.,    such  that 

\ a  \ +  Eajk  uj  (2) 

does  not  involve  applying  loads  to  the  added  constraints.    Thus,  for  k>p, 

P 
+  cX/i 


Define  the  filtering  matrix 


P 


0=1  J  J 


w  hich  is  such  that  ,  n    ,  ^  ,  * 

f  0    if  k£p 


Thus,  for  k>p, 


i^  if  k>p 


(3) 


F  ■  I  -  Su.uTm  W 


(5) 


CMPi^  -  CM^  =  u^vk  -  o^\  +  <  JaJkUj  (6a) 

J-l 
while  for  kip, 

CMF^  =  0  (6b) 

Equations  6  may  be  combined  to  give 

UA  +  B  »  CMFU  (7) 

where 

A  -  dlagCO  0  ...   0  uf  ^  oT22   ...  o^2]  (8) 
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and  B  is  a  matrix  the  first  p  columns  of  which  are  zero  and  the  remaining 


columns  of  which  are  linear  combinations  of  uu  Uj , 
the  matrix  \flB  =  U  MB,  the  kth  column  of  which  is 


u  .  Now  consider 


(9) 


and,  because  of  equation  lb,  this  can  have  nonzero  elements  only  in  the 
first  p  positions.  Thus  A  +  U-1B  and  (A  +  U"JB)n  are  of  the  forms 


l?2iDJ 


loj  ~rP~~  (10a>b) 

respectively,  where  Oj  is  a  p*p  matrix  of  zeros,  02  is  a  (N-p)*p  matrix 

of  zeros,  E  is  the  nonzero  part  of  U"lB,  and  D  is  the  nonzero  part  of  A. 
Finally,  from  the  equality 

(A  +  U"lB)n  =»  (U^CMFU)11  (11) 

by  taking  traces,  we  get 


N 

53"C2n    =  tr[(A+  tT'B)n]  -  tr[(ir»CMFU)n]  =  tr[(CMF)n]  =  trCcf)         (12) 
k=p+l 

where  Q  is  defined  as  the  matrix  triple  product  CMP.     Note,  however,  for 
computational  convenience,  that 

tr(cf )  =  tr[(CMF)n]  =  tr[(MFC)n]  =  tr[(FCM)n]  (13) 

We  have  used  the  fact  that  if  A  and  B  are  matrices  conformable  in  either  order, 
then  tr(AB)  =  tr(BA). 

V/e  will  refer  to  equation  12  as  the  nth  DM  (Dunkerley-Mikhlin) 
evaluation  and  to  the  corresponding  relation 
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vi^[tr(Qn)]"1/2n  (l4) 

as  the  nth  DM  approximation.     We  will  repeatedly  us«  the  symbol  £  to  mean 
"greater  than  but  approximately  equal  to,"  and  a  symbol  indicating  the 
opposite  order. 

If  M  is  diagonal,  if  there  are  no  rigid  body  modes  (i.e.,  if  p=0), 
and  if  we  take  n=*l,  equation  14  becomes        M 

a)j2  »    tr(Q)  =  tr(CM)  =  21\<^  (15) 

This  is  what  is  generally  referred  to  as  Dunkerley's  formula.     It  was 
(obliquely)  stated  In  189^  by  S.  Dunkerley  who  regarded  it  as  an  empirical 
representation  of  calculations  he  had  made  of  shaft  and  rotor  frequencies. 
Simple  extensions  to  continuous  mass  distributions  are  obvi6us  and 
have  been  discussed  by  many  writers,  among  the  first  of  whom  were  Temple 
and  Bickley,  Reference  6.       However,  the  theoretical  basis  for  the  contin- 
uous case  seems  to  have  been  first  established  by  Mikhlin,  who,  however, 
did  not  consider  problems  of  mechanical  vibration  and  who  seems  to  have  been 
unaware  of  Dunkerley's  work.     In  the  continuous  case,  matrix  multiplication 
is  replaced  by  integration  and  the  compliance  matrix  C  is  replaced  by  a 
symmetrical  function  of  two  variables,  z(x,y)  which  gives  deflection  at  x 

due  to  unit  loading  at  y.     This  notation  is  appropriate  to  a  one-dimensional 

herein 
field;  however,  examples  given  in  Roformnoo  1  illustrate  cases  involving 

two  dimensional  fields.     The  mass  matrix  M  which  heretofore  need  not  have 

been  a  diagonal  matrix,  now  becomes  the  equivalent  of  an  infinite  diagonal 

matrix,  namely,  a  function  m(x)  specifying  mass  per  unit  length.     The  matrix 

Q  is  replaced  by  the  function  T 

q(x,y)  «  m(y)[z(x,y)  -  2i^(y)    |     ii_(0  m(0  z(x,0  d£]  (16) 

k=l  K        ^0    K 

-6- 


where  the  functions  u.(x),  k=l,2,...,p,  describe  the  known  orthonormal 


rigid  body  modes,  satisfying 

L 


( 


m(C)  u±U)  UjCO  d£    =  &  (17) 


0 
(Kronecker  delta).     In  equation  16,  z(x,y)  is  the  compliance  of  any 

system  which  is  like  the  given  system  but  with  additional  constraints 

so  as  to  eliminate  the  rigid  body  modes.  The  function 

z(x,y)  is  sometimes  called  the  Green's  function. 

The  continuous  analog  of  raising  the  matrix  Q  to  the  nth  power 
is  the  formation  of  the  nth  iterated   function 


%(*,y) 


L 


q^CxjO  q(C,y)   d£;       qx(x,y)  =  q(x,y)  (lg) 


0 
and  the  operation  corresponding  to  taking  the  trace  of  Q    is 


y      -2n 

k=p+l 


q^Cx.x)  dx  (19) 

0 


If  there  is  both  lumped  and  distributed  mass,  we  may  write 

r 

m(x)  =  m(x)  +  2^  m    5(x-x.  )  (20) 

k=l    K  K 

where  6(x-x,  )  denotes  the  Dirac  "function,"  m(x)  denotes  a  continuous  mass 
distribution,  and  the  sum  represents  r  distinct  point  masses  m.  located 

aX>    X  a  Xi_,  K~J.,<i  |  .  .  .  ,r» 

Illustration 

We  consider  a  finite  element  model  of  a  segment  of  an  Euler-Bernoulli 

T  i 

beam.  The  displacement  vector  is  u  *  [Xj  x2  x3  xH]  where  xl  is  the  (up- 
ward lateral)  displacement  at  the  left  end,  x3  is  the  similar  displacement 
at  the  right  end,  x2/L  is  the  slope  at  the  left,  and  x^/L  is  the  slope  at 
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the  right.  The  length  of  the  segment  is  L.  The  FEM  consistent  stiffness 
and  mass  matrices  are 


K  =  (2EI/L3) 


'6     3-6     3' 

3    2-31 

-6  -3    6  -3 

L  3    1-3    2 


N  -  (m/420) 


156  22  54  -13 
22  4  13  -3 
54  13  156  -22 

__13  -3  -22       4 


First  we  exhibit  the  classical  solution.     The  modal  matrix  (the  first  two 
columns  of  which  are  not  uniquely  determined)  is 

U  -   Tl      a      b        cl/vfiT,     where  a2  ■  3,  b2  ■  5,  and  c2  =  7. 

0  -2a  -6b  -12c 

1  -a      b      -c 
0  -2a    6b  -12c 

and  the  spectral  matrix  is 

fl2  -  diag[u>2  w2  w2  u2]  *  (840EI/mL3)  diag[0     0     6/7    10] 

2       2        3       «♦ 

There  is  no  difficulty  in  verifying  equations  la  and  lb. 

Next,  we  employ  the  EM  method  and  verify  the  correctness  of 
equation  12.     Vfe  constrain  the  system  by  rigidly  fixing  the  left  end 
so  as  to  form  a  cantilever.     We  easily  determine 


C  -  (LV6EI) 


0 

0 

0 

o" 

0 

0 

0 

0 

0 

0 

2 

3 

0 

0 

3 

6 

CM 


(mL3/2520EI) 

"00   0   0 

0  0   0   0 

69  17  246  -32 

.84  21  336  -42. 

and  we  also  calculate 

P  - 

-3  -4   3  1" 

36  33  36  3 

3-1-3  4 

.36  3  -36  33. 

730 

Next,  we  obtain 

q  »  CMF  =  (mL3/25200EI) 

"0   0  0  0" 

0   0  0  0 

-3  -19  3  16  i 

0  -35  0  35 
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Finally  we  calculate 


tr(Q)  =  38mL3/25200EI  =  w3     +  M^ 


tr(Q2)  =  123M(mLV25200EI)2  «  u>:     +  a> 


-<♦ 


3  "♦ 

&  -6 


tr(Q3)  =  42902  (mLV25200EI)3  =  w3     +  u>„ 
and  so  on.     All  of  these  calculations  check  out  correctly. 

One  may  observe  that  the  two  diagonal  elements  of  Q  give  directly 


U). 


=  35mL3/25200EI,     u>~  =  3mL3/25200EI 


but  this  is  only  fortuitous  in  the  present  illustration  as  may  be  seen  in 
examples  lb  and  lc  which  appear  later. 

Now  we  illustrate  the  equations  employed  in  the  derivation  of 
equation  12.     First,  it  is  easy  to  verify  that 

FU 
Next  we  form 

CK  = 


[0     0    u     u  ] 

3       <♦ 


'  0     0 

0 

0" 

0     0 

0 

0 

-1  -1 

1 

0 

0  -1 

0 

1 

and  obtain 


[v3  vj  =  CK[u3  uj 


0 

o" 

0 

0 

6b 

10c 

12b 

0 

m 

/rfn 


which  obviously  satisfy  the  constraints  which  have  been  introduced  at  the 
left  end.     We  can  also  calculate  the  coefficients  a.l  3  =  2b,  a23  =  -3b/a, 
allf  ■  5c,  and  &2k  =  -6c/a,  so  that,  indeed 

v    =  u,  +  a,  ,u    +  a„u  •     v    =  u    +  a    u    +  a_   u 

3  3  131  232'  <♦  •♦ 


!<♦      1 


2>*     2 


We  can  also  calculate  and  verify  that 


u>2CMu,  =*  (6/7)(840EI/mL3)(bmL3/2520EI)[0  0  21  42]T/^n  =  v3 


with  a  similar  calculation  relating  u^  and  v^. 
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Next  we  calculate 

CMFU  »  QU  »   (•mtV^OEI) 


0 

0       0 

0 

0 

0       0 

0 

0 

0     7b 

c 

p 

0  14b 

0 

A  =  (mL3/25200EI)  diag[0  0  35  3] 


UA  -  (yBL3/25200EI) 


0  0      35b      3c 

0  0  -210b  -36c 

0  0      35b    -3c 

.0  0     210b  -36c 


B  -  CMFU  -  UA  -   (*4rlL3/25200EI) 


0  0  -35b  -3c 

0  0  210b  36c 

0  0  175b  33c 

0  0  210b  36c 


(mL,/25200EI)[0     0     35(a13Ul+a2  3u2)     3(3,^^)] 


Noting  that  U"1  -  (MU)   ,  we  finally  calculate 


U"XQU 


(mL3/25200EI) 


so  that 


D  =  (mL3/2520QEI)f35     6 
1  0     3 


0  0      70b    15c 

0  0  -3.5ab  -6ac 

0  0        35        0 

0  0  0        3  . 


E  =•   (mL3/25200EI) 


A  +  U~'B 


-ii 


70b      15c " 
-3.5ab    -6ac 


Vfe  also  verify  that 

tr(U"lQU)  -  38mL3/25200EI,     tr[(U_1QU)2]  -  123M(mL3/25200EI)2 ,  etc. 
Although  it  is  worth  remarking  that  neither  K  nor  C  possesses  an 
inverse,  this  is  to  be  expected  in  cases  for  which  there  are  rigid  body 
modes j  and  the  modified  system  involves  perfectly  rigid  constraints.     It 
is  of  some  interest  to  note  the  form  of  the  product  CK. 
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Applications 

Example  1,  part  a.  A  uniform,  massless  cantilever  beam  has  concentrated 
masses  m,  9m,  and  4m  at  distances  9L,  21L,  and  27L,  respectively,  from 
the  fixed  end.  It  is  desired  to  estimate  the  lowest  frequency  of  harmonic 
oscillation.  Taking  the  lateral  deflections  of  the  point  masses,  namely 
x,,  x2,  and  x^,  as  the  elements  of  a  vector  u,  we  easily  determine 
M  «  m  diagCl  9  *»];■ 


C  »  (9L3/EI) 


27  81  108 
81  343  ^90 
108  ^90  729 


Q  =  CM  =  (9mL3/EI) 


'  27  729  ^32 

81  3087  I960 

108  4410  2916 


H 


Prom  equation  12,  taking  n=l,  we  get 

V  V 

wl£(EI/54270mL3)  2=4.2926«10~3(EI/mL3)  2 

Taking  n=»2,  we  get 

u>,£  (E2I2/2878089084m2L6)^«  4.3174l6-10"3(EI/mL3) 
A  Rayleigh  approximation,  details  of  which  are  not  given  here,  gives 
the  numerical  coefficient  for  an  upper  bound  as  4. 317542* 10" 3,  so  that 
the  first  approximation  is  in  error  by  less  than  0,6%   and  the  second 
is  in  error  by  less  than  0.003$. 

Example  1»  part  b.  Now  suppose  that  the  fixed  support  at  the  left  of 
this  beam  is  replaced  by  a  frictionless  pivot.  The  frequency  cox  =  0 
corresponds  to  the  rigid  body  mode  u,  =  [3  7  9]T/(774)'/2  .  we  obtain 


F 


765  -189  -108 
-21  333  -252 
-27  -567    450 


/774;     Q  =   CMF  =  (9mL3/86EI) 


-702  -810  864 
-6198  -10962  10952 
-9858  -22806  21024 


and  we  have  no  difficulty  in  obtaining 

u3a^0.03195(EVmL3)^2,    wxS  0.033494 (EI/mL3)^ 
for  the  first  and  second  DM  approximations,  respectively.     A  Rayleigh 
approximation  gives  the  numerical  coefficient  as  0.033586  so  that  the 
errors  are  less  than  5%  and  0.3%  respectively. 


_n_ 


Example  1,  part  c.     Now  suppose  that  the  beam  is  completely  free  at  the  left 
so  that  there  are  two  rigid  body  modes,  involving  translation  and  rotation. 
We  use  the  same  matrix  C  as  above  although  a  simpler  one  could  be  constructed 

since  the  leftmost  segment  of  length  9L  obviously  does  not  enter  the  present 

T  \A  T 

problem.     Ws  take  ul  -  [ill]  /(14)  z  and  u^  as  a  combination  of  [1  1  1] 

T 
and  [3  7  9]    proportioned  for  orthonormality.     Thus  we  find  Uj  ■  ; 

[15  1  -6]T/(378)/2.     We  calculate 


42  -126     84 

-14      42  -28 

21    -63    42 


/126;     Q  -  (96mL3/EI) 


"0 

0  0' 

1 

-3  2 

L3 

-9  6J 

and  we  get  tr(Q)  «  288mL3/EI,  tr(Q  )  -  82944 (mL3/EI)2  both  of  which  give  the 
exactly  correct  result 

u>3  -  (EL/288mL3)^ 
since  there  are  only  three  frequencies  the  smallest  two  of  which  are  known 
to  be  zero. 

Example  2.     Now  we  consider  a  cantilever  beam  in  the  form  of  a  truncated 
right  circular  cone  having  base  radius  b,  tip  radius  a,  and  length  L,  and 
having  uniform  physical  properties.     Assuming  an  Euler-Bernoulli  model, 
we  obtain  the  bending  compliance 

z(x,y)  -  {(£-n)[(l-an)"2-  (l+2an)]  +  2a2n7(l-an)>/6Ba2 
where 

B  =  irEbV4  =  (EI)root;     a  3  (b-a)/bL;     £  =*  Max{x,y};  n  =  Min{x,y} 
The  mass  per  unit  length  is 

m(x)  =  mo(l-ax)2;      m©  ■  irb2Y 
mo  being  the  value  pertaining  to  the  root  end  and  y  denoting  mass  density. 
Performing  the  indicated  evaluations  gives 

W*;2  £      m(x)  z(x,x)  dx    =  moL"(l+4q)/60B 
for  the  first  estimate  and,  noting  the  symmetry  of  z(x,y), 
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in  it 


m(x)  m(y)   z  (x,y)  dx  dy  = 

=  nv,2L8[3(q,,-l)/8  -  qkloge(q)/4  +  3p2  -  5p3  +  llpV4 

-  3p5/5  -  p6/20  -  p7/70  +  9p8/80  -  p9/5  +  2p» °/25]/l8B2p8 
where     q  =  a/b,  p  =  1-q.     Although  both  of  these  evaluations  give  results  for  u^ 
which  are  smaller  than  the  correct  value,  the  second  gives  a  result  which 
has  a  maximum  error  of  0.8%  (for  q  =  0)   for  0  =  q  =  1.  The  first,  and 
much  simpler,  evaluation  is  much  less  accurate,  but  it  serves  as  a  guide  to 
the  following  approximation 

w15-  8.72[(l-0.0l6q)/(l+5.053q)]/2(B/moL")/2 
which  has  a  maximum  error  of  \A%  for  q  =  0.1.     The  "exact"  results,  for 
comparison,  were  calculated  from  formulas  given  by  Conway,  Becker,  and 
Dubil,  Reference  2. 

Next,  we  consider  some  cases  for  which  the  mass  distribution  is 
of  the  form  described  by  equation  20.  The  first  and  second  DM  evaluations, 
respectively,  take  the  forms 


Z-)u"2=        m(x)  z(x,x)  dx    +    /-jm.    z(x,  ,x,  ) 
i=l  1       *o  k»l  *        K    K 


Z^w""=  m(x)  m(y)   z2(x,y)  dy  dx    +  2  ^^ m.        m(x)   z  (x^)  dx 

1*1         Vo  k=1       0 

+  2-,  2-im,  m.  z  (x.  ,x.) 
i=l  j=l         J  J 

The  first  term  in  these  evaluations  is  what  we  had  previously  for  continuous 
distributions  and  the  last  term  may  be  identified  with  the  trace  of  the 
appropriate  matrix. 
Example  3.     Suppose  that  the  cantilever  beam  of  example  la  itself  has  a 
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mass  4m  distributed  uniformly  along  its  length.   (Note  that  the  upper 
limit  of  integration  is  27L. )     We  have  m(x)  =  4m/27L  and  z(x,y)  ■ 
(3xyw-w  )/6EI,  where  w  ■  Min{x,yJ.       We  calculate 

w,  £  [(W27L)(27L)Vl2EI  +  54270mL7EI]"/*  =  0.004054(EI/mL3) ^ 
for  the  first  DM  approximation.     Note  that  we  have  seen  the  number 
54270  in  example  la.     The  second  DM  approximation  is  more  labor;  it 
gives 
Uj  £  [ll«3,7/35  +  4759166988/7  +  2878089084]"^ (EI/mL3)1'* 
=  (81O.554923033/35)"  ^(EI/mL3)*2  =  0. 004  08289  (EI/mL3)^2 
An  "exact"  evaluation  by  transfer  matrix  procedure  gives  the  numerical 
coefficient  0.00408305,  so  that  th6  errors  for  the  first  and  second 
approximations  are  0.7%  and  0.004$  respectively. 

Example  4.     Consider  the  case  of  a  uniform  cantilever  beam  of  length  L 
and  mass  pm  having  a  concentrated  mass  (l-p)m  at  its  tip.     The  compliance 
function  z(x,y)  is  given  in  example  3    and  the  mass  function  is  m(x)  a 
pm/L  +  (1-p)  m  6(x-L).       The  first  and  second  DM  approximations  give 
u)^  [(3EI/mL3)/(l-3p/4)]^;  {(a):-1.5?  at  p  -  1} 
w^  [5040(EI/mL3)V(560-856p+329p2)]/S     {(b):-0.017*  at  p  -  1} 
( In  these  and  other  estimates  given  in  this  example  4 ,  the  estimate  is 
distinguished  by  a  lower  case  letter  identification,  followed  by  the  maximum 
percent  error  and  the  value  of  p  for  which  it  obtains.) 

This  example  affords  an  opportunity  to  discuss  other  types  of  simple 
but  accurate  approximations.     A  Rayleigh  approximation  based  upon  the 
deflection  function  y  »  3Lx2-x3  gives 

Wj  4  [(3EI/mL3)/(l-107p/l40)]^;     {(c):+1.5*  at  p  -  1} 
This  is  an  upper  estimate  while  (a)  is  a  lower  estimate.     An  obvious  average 
is 

u.i  »  [(3EI/mL3)/(l-53p/70)]'/z;         {(d):-0.392  at  p  -  0.83  j 
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Approximation  (d)  gives  equal  weighting  to  (a)  and  (c).     This  happens 
to  be  approximately  optimal  in  this  case  but  this  is  only  fortuitous.     The 
coefficient  of  p  should  lie  somewhere  between  105/1^0  and  107/1^0.     Both  are 
exact  for  p  =  0.     We  can  choose  the  coefficient  so  as  to  make  the  result 
exact  for  another  value  of   p ,  say  p  =  1.     This  is  a  uniform  cantilever  for 
which 

u)l  »  (1. 87510^0687... )2(EI/mL3)K2 
Thus  we  obtain 

w1  »  [(3EI/mL3)/(l-0.757328p)y2;       {(e):-. 3662  at  p  -  0.81} 
Exactly  the  same  result  is  obtained  by  replacing  the  first  term  in  the 
first  DM  evaluation  by  the  known  exact  value  appropriate  for  a  uniform 
cantilever  with  no  added  mass  at  the  tip.     This  suggests  replacing  the 
coefficient  11/1680  in  the  first  term  of  the  second  DM  evaluation  by 
(1.8751...)"8  and  we  get 

a),  js   [630(EI/mL3)2/(70-107p+4l.l222805p2)]lyS     {(f):-0.006#  at  p  =  0.86} 
This  accuracy  is  almost  incredible.     It  goes  without  saying  that  the 
accuracy  of  these  approximations  is  far  greater  than  that  of  the  physical 
theory  to  which  they  pertain. 

Example  5.    we  next  take  up      an  example  which  combines  the  features  of 
mixed  (i.e.,  lumped  and  distributed)  mass  distribution  and  rigid  body 
motion.     We  consider  the  axial  motion  of  a  uniform  elastic  bar  of  mass  m, 
length  L,  and  axial  stiffness  EA,  which  has  a  concentrated  mass  2m  at  the 
left  end  and  a  concentrated  mass  3m  at  the  right  end.     This  dumb-bell 
shaped  object  is  not  constrained  or  tied  down  at  any  point.     The  mass 
function  is 

m(x)  ■  m/L  +  2m6(x-0)  +  3n6(x-L) 

We  consider  the  compliance  of  a  system  fixed  at  the  left  end. 

z(x,y)  =  Min(x,y)/AE 
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and  we  obtain 

q(x,y)  =  [m(y)/AE][Min(x,y)  -  (8Lx-x2)/12L] 
Thus  the  first  DM  evaluation  gives 

a),  >  (9AE/13mL)^  =  0.832(AE/mL)^ 
and  the  second  gives 

u>   >(l^A2E2/242m2L2)^  =  0.87829(AE/mL/2 
The  exactly  correct  coefficient  is  the  smallest  positive  root  of  the 
equation  tan  x  =  5x/(ox  -1),  which  is  approximately  0.87935.     Thus  the 
errors  are  5»^%  and  0.12$  respectively. 


Example   6     We  determine  the  fundamental 

harmonic  frequency  of  a  uniform  hinged- 

, „   ,  *_  _«_  T       ^  *  Sketch  of  uniform 

guided  beam  of  length  L.     Tne  mass  func-  hinged-guided  beam 

tion  is  m(x)  ■  mo  ■  const.     The  compliance 

function  is 

z(x,y)  *  z(y,x)  =  [3xy(2L-y)  -  x3  +  <x-y>3]/6EI 


'(x-y)3  if  x*y 
0  if  x£y 


<x-y>3  - 


1 


The  DME1,  given  by  equation  19  with  n  »  1,  is 

/•L 


u> 


f  £   S^2  =   [  moZ(x,x)dx  =  moLV6EI;  wx£  2.^95(EI/m0L-/* 


The  DME2  is  obtained  from  equation  19  with  n  *  2,  viz. 

^  £  Z<*?  =     f  f  q(x,5)q(5,x)dSdx  = 

i=l  1         >oh 


»  2I  f  m(x)m(Oz2(x,Od£dx 


>(  (  m(x)ra( 
'0*0 

because  of  the  symmetry  of  z(x,y).  Tne  important  consequence  of  this 

is  that  we  can  take 

z2(x,0  =  [3x£(2I^x)-53]2/36E2I2 

and  not  have  to  concern  ourselves  with  the  alternate  form  applicable 

for  £  >  x.     Continuing  the  evaluation  we  arrive  at  (17/630) (moL1* /EI)2, 
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from  which  we  obtain 

WjSs  2.W3(EI/m0l/,)Vk 

The  exact  result  is 

wx=  (Tr2A)(EI/m0L,#),/Z=  2.M67*l(EI/m0L,|)" 

so  that  the  excellence  of  the  approximations  is  evident. 

Having  the  excellent  approximation  to  uix   given  by  the  DMA2,  we 


can  get  a  reasonable  approximation  to  the  second  frequency,  viz. 

*DMA1  "  ^l    )EMA2 


ufcz   (w72)mfll  -  C(^)mA^y'=  ^6  "  (17/630 ),/e]moL-/EI 


giving 


o),«20.4(EI/moLH)Va 


2 

The  correct  coefficient  is  (3V2)2  sr  22.2. 

This  exanple  permits  us  also  to  exhibit  another  possibly  useful 
feature  of  the  method.  In  this  case  it  is  known  that 

\  -  [(k-iM^EZ/taoL*)*,  k=l,2,..., 

-2      -4 

so  that  our  precise  evaluations  of  Za>.  and  Eu),  lead  to 

ao  co 

22    n"4  =  irV96;    2j   n"8  =  17*7l6l280 
n=l,3,5,...  n=l,3,5,... 

These  are  well  known  results;  cf.  Jolley  (5).  However,  for  other  cases 

the  corresponding  sums  may  be  unknown.  For  this  reason  we  have  looked 

into  all  possible  cases  of  a  uniform  beam  having  end  conditions  of  the 

following  types:  fixed,  hinged,  guided  (see  right  end  of  the  figure),  and 

free.  Letting 

n*  =  u>2m0L,*/EI 

the  results  may  be  exhibited  in  the  form  shown  in  Table  1. 

In  Table  1  the  physical  case  described  by  the  specification  of  end 

conditions  corresponds  to  the  frequency  equation  indicated.  The  sum  of 

reciprocal  second  and  fourth  powers  of  the  frequencies  (except  for  factors 
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Table  1   Some  results  for  uniform  beams 


CASE 

EQUATION 

END  CONDITIONS 

r  -2 

*■* 

Hi 

1 

cosn  s  sechn 

clamped-cl amped ( * ) 
free-free 

1/420 

71/17463600 

4.7300407449 

2 

tann  *  tanhn 

clamped-hinged(  *) 
hinged- free 

1/210 

13/727650 

3.9266023120 

3 

tann  =  -tanhn 

clamped-guided ( * ) 
guided-free 

1/30 

29/2835 

2.3650203724 

4 

cosn  s  -sechri 

clamped-free 

1/12 

11/1680 

1.8751040687 

5 

slnn  m  0 

hinged-hinged ( * ) 
guided-guided ( * ) 

1/90 

1/9450 

3.1415926536 

6 

cosn  a  0 

hinged-guided 

1/6 

17/630 

1.5707963268 

*  The  root,  n  3  0,  is  to  be  excluded  for  cases  marked  with  an  asterisk 


of  El/m ,Lk   and  E2I2/mo2L8,  respecively)  are  shown.  The  smallest  root 
of  the  frequency  equation  is  also  given.  This  information  may  be  useful 
in  augmenting  information  contained  in  work  by  Young  et  al.  (16)  and  by 
Gorman  (4).  It  may  be  remarked  that  neither  of  these  references  treats 
the  end  condition  ^  =  ^-   =  0  which  we  call   "guided;"  this  condi- 
tion is  certainly  not  technically  important  but  its  inclusion  seems  to 
be  indicated  on  the  grounds  of  completeness. 

In  the  table,  an  asterisk  indicates  that  the  obvious  root,  n  -  0, 
Of  the  frequency  equation  is  to  be  excluded.  The  sums  include  only  non- 
zero frequencies,  and  nl   Indicates  the  smallest  positive  root.  We 
evidently  have  such  results  as 


d) 


ZCn?5]'*-  71/17463600 
i-i 


where  n.^  denotes  the  ith  positive  root  of  equation  j.  From  this 
evaluation  we  could  obtain 

n ^4. 7191 
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Exanple  7     In  order  to  illustrate  a  case  for  which  the  field  of  integration 

is  two  dimensional,  we  consider  a  uniform  circular  elastic  plate  which  is 

clamped  at  its  outer  radius  r  =  a.       The  influence  relationship,  z(P,Q), 

is  given  by  Tlmoshenko  (Ik)  in  the  form  of  a  series,  but  for  our  purposes 

all  that  is  needed  is  the  "self"  influence,  i.e.,  the  deflection  under  the 

unit  lateral  load 

z(P,P)  =  (a2-r2)2/l6TTDa2 

which  is  also  given  in  (Ik).     Thus,  with  mo  denoting  mass  per  unit  area  of 

plate,  we  find  ^ 

a>72  £  zLcjT2    =(mo/l67TDa2)f        (a2-r2)2rdrd8  =  i^aV^D 
~   i-1  jo  'n 


whence  ,. 

w,^  6.9(D/moa*)* 


This  is  actually  a  poor  estimate,  the  correct  coefficient  being 
approximately  10.22  rather  than  6.9.  The  reason  lies  essentially  in  the 
fact  that  the  sum  includes  a  double  infinity  of  terms  corresponding  to 
both  radial  and  circumferential  nodal  lines  and  the  frequencies  are  not 
well  separated.  A  better  approximation  may  be  obtained  with  great  labor 
by  considering  the  EME2,  or,  much  more  easily,  by  excluding  all  modes 
having  radial  nodal  lines.  This  is  accomplished  by  considering  a  unit 
load  uniformly  distributed  on  the  circumference  of  a  circle  of  radius  p. 
The  deflection  is  (1*1) 

z(r,p)  =  [(r2+p2)loge(p/a)  +  (l+p2/a2)(a2-r2)/2  -  (p2-r2)]/8TrD 
for  r  £  p,  and 

z(r,p)  =  [Cr2+p2)loge(r/a)  +  (l+p2/a2)(a2-r2)/2]/&VD 
for  r  £  p.       The  mass  function  is 

m(r)  =  27rrmo 
and  without  difficulty  we  obtain  the  DME1  and  the  DME2 
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eo  oo 

2w72    -  moaV96D;     £<*?*  s  (17/184320)  (moaVD)* 
i»l  x  1=1  x 


from  which  we  get  the  respective  estimates 

uijs:  9.8(0/1^")*      u)1s:10.2042(D/moa,,)Va- 

It  is  difficult  to  assess  the  accuracy  of  the  latter  result. 

Ihe  correct  coefficient  is  the  square  of  the  smallest  positive  root 

of  the  equation 

J0(x)  IjCx)  +  Jj(x)  I0(x)  =  0 

which  is  given  by  Rayleigh  (9)  as  (3.20)  »  10.24.  However,  an  upper 

bound  given  by  Timoshenko  (13)  is 

[19.2(51-2(519)'  ]*  «  10.217 

so  that  our  result  is  definitely  in  error  by  no  more  than  0.13$. 

Example  8   Consider  a  uniformly  tensioned 
uniform  circular  elastic  membrane  of  radius 
a.  The  influence  function  is 

z(P,Q)  =»  z(r,e;p,<|>)  = 

*  ci/ihrT)  ioge{a-^p2:2a2rpco;(9'f)j 

^  a2[r2+p2-2rpcos(M>)] 
which  is  a  rearrangement  of  a  form  given  by 

Rektorys  (10)  with  the  addition  of  the  membrane 

tension  T  in  the  denominator.     Ihe  DME1  gives 


Notation  for  analysis 
of  circular  membrane. 


i=l 


-2 


'z(P, 


P)dS  =    ~ 


Thus,  the  DME1  is  useless  except  that  it  demonstrates  the  divergence  of 
the  double  sum 

m=0  n=l  »n 
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where  j      .  is  the  nth  positive  zero  of  J  (x).     We  have  been  unable  to 
m,n  m 

evaluate  the  integral  which  gives  the  DME2.     For  a  successful  numerical 
evaluation  we  do  as  was  done  in  the  preceding  exanple,  namely  confine 
attention  to  axisynmetric  modes.     The  deflection  at  radius  x  of  such  a 
circular  membrane  loaded  by  unit  force  uniformly  distributed  on  a 
circle  of  radius  y  is 

z(x,y)  -  (1/2ttT)  loge[a/Max(x,y)] 
and  the  mass  function  is 

m(x)  =  2TTmoX 
There  is  no  difficulty  in  obtaining  the  DMEt  and  the  DME2,  viz. 

co         -  oo 

2>7   =  moaVMT      23  w74  =  raJaV32T2 
i«l  1  i=l  ± 

As  in  the  evaluation  in  the  preceding  example,  these  sums  include 
only  frequencies  for  axisymmetric  modes.  The  exact  results  in  the 
present  case  are 

Since  j0  ,»  2.4048,  the  DMA1,  which  gives  the  coefficient  2,  is  in  error 
by  17%   and  the  DMA2,  which  gives  the  coefficient  (32)/if  is  in  error  by 
only  1,135. 

Incidentally,  we  have  here  established  the  interesting  and  possibly 
useful  results  . 


2j-2.  =  1/4;   2jA  =  1/32 
i»l  u»1         i=l  u>1 


Example  9  We  now  consider  the  lateral  vibrations  of  a  simply  supported 
rectangular  elastic  plate  having  dimensions  a  and  b.  V/e  assume  that  the 
ratio  p  ■  b/a  does  not  exceed  unity.  The  compliance  is  (14) 

00  oo 

z(x,y;£,n)  =  (V^abD)^  ^[sindmrx/a)  sin(rmr^/a)  sin(rory/b)  sin(mrn/b)] 

m=l  n=«l  [(m/a)2  +  (n/b)2]2 
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Thus,  the  DMEL  is 

2>![2=  (moa'A))^  £  (™2*2  +  nVp2)'2 
lal  m=l  n-1 

»  (m0aV1*TT,,p,|Dy    [n27r2p2csch2(mrp)  +  mrpcoth(mrp)  -  21/n* 
n=*l 
which  converges  fairly  rapidly.     The  WE2  gives 

qq  00        oo 

2w^  =  (iHbaVir*D)aV  ^T  (m2  +  n2p2)'* 

m=l  n=l 

CO 

*  (mgaV96peTT8D2)/  [2nf>TTltp't(csch'>nTTp+2csch2nTrpcothnTTp) 
n=I 

+  12n3TT3p3csch2mrpcothmTp 

15mrp<n7rpcsch2mTp+cothmrp-A  8]/n 8 

and  this  series  converges  quite  rapidly.  Thus,  for  the  particular  case 
a  ■  b,  by  taking  only  one  term  of  the  series  we  get 

(Sf  «  6.7556«10"*6m2a8/D2;a)l«19.6l5(D/moa,>)/2 
The  exactly  correct  coefficient  is  2tt2  »  19.739  so  that  our  estimate  is 
quite  good. 

For  this  particular  problem,  a  simple  formula  (17)  gives  all  of  the 
double  infinity  of  natural  frequencies,  viz. 

wm  n  "  1T^m2/a2  +  n2/b2)(D/m0)/i 
and  the  correctness  of  the  DI^El  and  the  DP-E2  can  be  verified  term  by  term. 
It  is  merely  fortuitous  that  the  exact  result  for  u;,  (=t*\i)  can  be  obtained 
in  this  case  by  taking  only  one  term  of  the  double  series  shown  above. 
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